JASR 6703 
TYPED/ 3/4/04 


ARTICLE IN PRESS 


No. of pages: 8 
DTD 4.3.1 / SPS-N 



ELSEVIER 


Available online at wvyfw.sciencedirect.com 

f 

• CIINCI //YlDIRICTe 




Advances in Space Research xxx (2004) xxx-xxx 


Advances in 

Space 

Research 

( a COSPAR publication) 


wyvw.elsevier.com/locate/asr 


2 An improved Green’s function for ion beam transport 

3 J. Tweed a ’*, J.W. Wilson b , R.K. Tripathi b 

4 * Department of Mathematics and Statistics, Old Dominion University, Hampton Blvd, Norfolk, VA 23529-0077, USA 

5 b NASA Langley Research Center, Hampton, VA 23681-2199, USA 

6 Received 19 October 2002; received in revised form 16 November 2003; accepted 23 November 2003 


7 Abstract 

8 Ion beam transport theory allows testing of material transmission properties in the laboratory environment generated by particle 

9 accelerators. This is a necessary step in materials development and evaluation for space use. The approximations used in solving the 

10 Boltzmann transport equation for the space setting are often not sufficient for laboratory work and those issues are the main 

1 1 emphasis of the present work. In consequence, an analytic solution of the linear Boltzmann equation is pursued in the form of a 

12 Green’s function allowing flexibility in application to a broad range of boundary value problems. It has been established that simple 

13 solutions can be found for the high charge and energy (HZE) by ignoring nuclear energy downshifts and dispersion. Such solutions 

14 were found to be supported by experimental evidence with HZE ion beams when multiple scattering was added. Lacking from the 

1 5 prior solutions were range and energy straggling and energy downshift with dispersion associated with nuclear events. Recently, we 

16 have found global solutions including these effects providing a broader class of HZE ion solutions. 

17 © 2004 COSPAR. Published by Elsevier Ltd. All rights reserved. 

18 Keywords: Radiation risk; Boltzman transport equation; Ion beam transport; An improved Green’s function 


19 1. Introduction 

20 In space radiation transport, the energy lost through 

21 atomic collisions is treated as averaged processes over 

22 the many events which occur over even relatively small 

23 dimensions of most materials and is referred to as the 

24 continuous slowing down approximation. It is reasoned 

25 that the few percent energy fluctuation in energy loss has 

26 little meaning for ions of broad energy spectra and es- 

27 pecially in comparison to the many nuclear events for 

28 which uncertainties are still relatively large. In contrast, 

29 the laboratory testing of potential shielding materials 

30 uses nearly monoenergetic ion beams in which the 

31 interpretation of the interaction with shield materials 

32 requires a detailed description of the interaction process 

33 for comparison to detector responses (Schimmerling 

34 et al., 1986). The development of a Green’s function 
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approach to ion transport facilitates the modeling of 35 
laboratory radiation environments and allows for the 36 
direct testing of transport approximations of material 37 
transmission properties. For a number of years, this 38 
approach has played a fundamental role in transport 39 
calculations for high-charge high-energy (HZE) ions 40 
and has been used to great effect by radiation investi- 41 
gators at the NASA, Langley Research Center. These 42 
earlier works have not, however, taken into account 43 
such effects as straggling or of the energy downshift with 44 
dispersion which occur whenever a nuclear event takes 45 
place. In addition to the validation of physical processes, 46 
a theoretical model of the role of straggling is essential 47 
to understanding of the radiobiology of ion beams as 48 
required in evaluation of astronaut risks which must be 49 
mini miz ed at least to within some regulated level (Shinn 50 
et al., 1999). The present development is in the context 51 
of an asymptotic expansion of the 3D Boltzmann 52 
equation, for which, the lowest order term is along the 53 
forward ray. Additional asymptotic terms are discussed 54 
in an earlier work (Wilson et al., 1991) and a related 55 
paper (Wilson et al., 2002a). 56 


0273-1 177/S30 © 2004 COSPAR. Published by Elsevier Ltd. All rights reserved. 
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57 2. The Boltzmann equation 


58 The specification of the interior environment of a 

59 spacecraft and evaluation of the effects on the astronaut 

60 is at the heart of the space radiation protection problem. 

61 For some time investigators at The NASA Langley 

62 Research Center have been developing techniques to 

63 address this problem and an in-depth presentation of 

64 their work is given by Wilson et al. (1991) although 

65 considerable progress has been made since that publi- 

66 cation (Cucinotta et al., 1998). The relevant transport 

67 equation is the linear Boltzmann equation. The lowest 

68 order asymptotic term is the straightahead approxima- 

69 tion. With the target secondary fragments neglected, 

70 Wilson et al. (1991), this equation takes the following 

71 form: 

S z <t>j(Z,E) = J2 [ <yAE,E‘)<i> k (z,E')&E l 

k>j J 

-a j(E)<Pj(Z,E), z>Z', ( 1 ) 

73 where <pj{z,E) is the flux of ions of type j moving along 

74 the z-axis at energy E in units of MeV/amu and oj{E) 

75 and o jk {E,E') are the media macroscopic cross-sections. 

76 The <jj k (E,E r ) represent all those processes by which 

77 type k particles moving in the z-direction with energy E 

78 produce a type j particle with energy E moving in the 

79 same direction. Note that there may be several reactions 

80 which produce a particular product, and the appropriate 

81 cross-sections for Eq. (1) are the inclusive ones. The 

82 total cross-section Oj(E) with the medium for each 

83 particle type of energy E may be expanded as 

°j(E) = o*{E)+of(E) + o r J (E), (2) 

85 where the first term refers to collision with atomic 

86 electrons, the second term is for elastic nuclear scatter- 

87 ing, and the third term describes nuclear reactions. The 

88 corresponding differential cross-section is given as 


a Jk (E,E') = Y^^nWjkSiE-E' + s n ) 


r {E + X jk -E'f} 


x exp 


A* 


( 3 ) 


90 where e„ are the atomic/molecular excitation energy 

91 levels and where the collision energy downshift X jk and 

92 corresponding energy width e jk are approximated from 

93 the known momentum distributions observed in heavy 

94 ion reactions and represented by a gaussian model. 

95 Many atomic collisions (~10 6 ) occur in a centimeter of 

96 ordinary matter, whereas ~10 3 nuclear coulomb elastic 

97 collisions occur per centimeter, while nuclear reactions 

98 are separated by a fraction to many centimeters 

99 depending on energy and particle type. This ordering 


allows flexibility in expanding solutions to the Boltz- 
mann equation as a sequence of physical perturbative 
approximations. 

We require to solve Eq. (1) subject to a boundary 
condition of the type </> y (z',£) = Fj(E). In the case of a 
unit source at the boundary, F}(£) takes the special form 

Fj(E) = S jk S(E - E'), (4) 

and the corresponding solution, which is called the 
Green’s function, is denoted by the symbol G jk {z,z!,E,E'). 
Once the Green’s function is known the solution for an 
arbitrary boundary condition Fj(E) is then given by 

4> j (z,E) = J2 [ G Jk (z,z!,E,E")F k (E")dE". (5) 

In the case of an accelerator beam, the boundary condi- 
tion consists of a narrow gaussian function in energy and 
is incorporated by addition to the straggling width on 
leaving the boundary. In the case of space radiations, the 
boundary condition is represented as a broad function of 
energy and direction for each ion type and is handled by 
ordinary numerical procedures. It should also be noted 
that Eq. (5) provides a basis for multiple layers of mate- 
rials by matching the solution at the boundary interface. 


3. Solution methods 


We rewrite Eq. (1) in operator notation by defining a 
vector array field function as 

<P = [4>j(z,E) ], (6) 

the drift operator 

D = [6J, (7) 


the interaction operator 


I = 3-<t = 


j o jk {E,E')AE! -[*,(£)], 


( 8 ) 


with the understanding that I has three parts associated 
with atomic, elastic, and reactive processes as given in 
Eqs. (2) and (3). Eq. (1) is then rewritten as 

£> • = / • <P = [/“ +r'+/ r ]-4>, (9) 

and one must look for solutions. In what follows, we 
will recall the solution of the atomic interactions by 
Payne (1969) and implemented by Wilson et al. (2002b). 
Effectively, we look at 

D<P = F'<P, (10) 


which must then be coupled to the remaining terms in 
Eq. (9). For analysis, it will be advantageous to make 
the following separations: 


[d - / at - r 1 + o r \ ■ $ = 


j v T jk (E,E')dE' 


<P = 3 r <P. 
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142 3.1. Atomic processes 


143 The lowest order approximation to the Boltzmann 

144 equation is given in terms of the atomic collision pro- 

145 cesses as 

D <P = / al <f> , (12) 

147 with the boundary condition 

*B = [ = (13) 

149 The solution, which incorporates energy straggling, 

150 takes the form 


V2rt k (z-E) P 2E k (z - E) 2 


152 where 

(E! k (z - E)) = R k ' [R k (E') - (z - z')], (15) 

154 where R k (E) is the usual range-energy relation and 

155 E k (z - E) is the rms deviation for incident &-type parti- 

156 cles of energy E 1 after a distance of penetration z — E 

157 (Wilson et al., 2002b). 

158 3.2. Elastic scattering processes 

1 59 The addition of elastic scattering processes is given by 

D <£=[/*' +/ el ]•*. (16) 

161 Since we have approximated the elastic scattering dis- 

162 tribution by 

<r«(E,E')=of(E')d Jk 6(E-E'), (17) 

164 we find that 

l'“]*«[0], (18) 

166 and thus 

D 4> »[/«]•#. (19) 

168 Elastic scattering does not appear in the first asymptotic 

169 term evaluated herein. The first correction will contain 

170 elastic scattering as a dominant term for the propagation 

171 of the surviving primary beam ions and in some bound- 

172 ary problems involving collimators elastic scattering will 

173 play a role for higher order terms. The elastic scattering 

174 propagator is a focus of current research and will couple 

175 with the present formalism. In the past, this coupling was 

176 in terms of acceptance functions and provided good 

177 agreement with neon ion beams (Shavers et al., 1993). 


[D - 7 a ‘ + <r r ] • tf> = J a T Jk {E, E') d£'J • 

= S r • <P. (20) 

In the present work, we approximate the fragment en- 181 
ergy distribution by 182 

gj i E exD f (-E + -v.-g)1 , 2n 

where X jk is the collision energy downshift (MeV/amu) 184 
and tj k is the interaction energy width (MeV/amu). Xj k is 185 
related to the momentum downshift (MeV/c) 186 


i22) 

via the equation 188 

(23) 

Aj(m + E) 

where A k is the projectile mass (amu), A j is the fragment 190 
mass (amu), E is the fragment energy (MeV/amu), m is 191 
the energy equivalent of a proton mass and 192 

p{E) = s/E 2 + 2mE. (24) 

is the fragment momentum (MeV/amu/c). The interac- 194 
tion energy width is similarly related to the momentum 195 
width «7 F (MeV/c) through the equation 196 

p{E)o F 

(25) 

where <r F is given as (Tripathi et al., 1994) 198 

1 / 45 25 \(A J {A k -Aj)\ 

We start with the solution of the equation 200 

\D — 7 at + <r r ] • G° = (0], (27) 

for a unit source at the boundary. Note that G° is di- 202 
agonal and takes the form 203 

J r c>\ _ ^r(E) 5jk 


(?Jz,E ,E,E') = 


Pj{E) \/2nSt(z — E) 


x exp - 


(E-jE'^z-E))) 1 

2s' k (z-E) 2 


where the nuclear attenuation is described by the func- 205 
tion 206 


Pk{E) = exp 


i-rs-]' 


178 3.3. Nuclear reactive processes 

179 Following the above analysis, we are left with 


and S k (E) is the change in E per unit path length per 208 
nucleon. Eq. (28) and the reactive integral operator are 209 
all that is required to develop the solution under the 210 


JASR 6703 
TYPED/ 3/4/04 


ARTICLE IN PRESS 


No. of pages: 8 
DTD 4.3.1 / SPS-N 


211 

212 

213 

214 

215 

216 

217 

218 


220 

221 

222 


224 

225 

226 

227 

228 

229 

230 

231 

232 

233 

234 

235 

236 

237 


240 

241 


243 

244 

245 

246 

247 

248 


4 


J Tweed et al. I Advances in Space Research xxx (2004) xxx-xxx 


straightahead approximation. The lateral spread of the 
beam is beyond the scope of the present development. 
So far all of the operators have had only diagonal ele- 
ments. Off-diagonal elements enter through the reactive 
regeneration terms <f Jk which appear on the right side of 
Eq. (20), The challenge is to further develop the solution 
of Eq. (20) and this will be accomplished as follows. The 
integral form of Eq. (20) can be written as 

#=(£> — / at + a r r 1 •< P B + j [D-I“ + or'-S' $dz l 

= (7° • <&„ + O • (7° • 3 1 ■ <b (30) 

where 4> B is the appropriate boundary condition. Eq. 
(30) is a Volterra integral equation and is easily solved in 
a Neumann series as 

<p = [& + QGPS r <? + Q(?-5 I QCpE r 
■ Q ■ (7° + ■■•]• #b 

= [G° + G 1 + G 2 -I ] • 4>b, (31) 

with the elements of the leading term given as Eq. (28). 
The above formalism lends the following interpretation 
of the solution. The operator G° propagates the particles 
with attenuation processes. The first term G° ■ 4> B 
propagates the ions at the boundary to the interior. 
S T ■ G° • is the production density of first generation 
secondaries at depth z\. These are propagated to the 
interior by G° • a* ■ G° ■ <f> B - Lastly, G 1 • 4> B = Q 
G° • a 1 ■ G° • 4> b represents the sum of all the first gen- 
eration secondaries being propagated from the interval 
[/, z] and so on. We have already identified the propa- 
gator G°. We now need to identify the remaining terms 
in the Neumann series and we begin noting that these 
are related via the recurrence formula 

G" +1 =[0-G°-,S r ]-G", w^0. (32) 


3.4. First collision term 

The second term in Eq. (31) is the first collision term 
G) k (z,z',E,E') = [Q • G° • S r • G\(z,z>,E,E) 

= JJ <?jM-zuE,E 2 ) 

* { / 

xG2 t (z,z',£ > £')dEiJdE 2 dzi. (33) 

The physical interpretation is that S r • G° is the volume 
source of ions from collisions at z\ of a unit ion source at 
z' of energy E' The ions present at z with energy E are the 
result of propagation from the all the ions through out 
the volume. The first task is to evaluate the volume 
source term 


G°yz 1 ,z',E 2 ,£') 




f ^k(Ei) 

J J(2n)e/ XP 

(Ei + A jk — E[) 2 
2e /* 


PfiE 1 ) 1 


\ (Ei — (E’fiz — z!))) 1 

PfiEi) fifins'fiz - 

-z') CXP 

2s> k {z-z!) 2 


d Ei. 


(34) 


Note that a sharp maximum occurs at = (E' k (z - z 1 )), 
E 2 = E\- k jk and the cross-sections and attenuation 
functions are siowiy varying functions of energy so that 
Eq. (34) can be accurately approximated as 


(E r • G°]^(zi ,z!,E 2 , E') 

Pk{E) 

P k {(E' k {zy -z'J}) 

1 


<*(», -m 


x exp 


{- 


( E 2 + lj k ~ {E’ k (z - z )) f 

2 l s 'k( z ~ + £ y*l 

The next step is to construct the term 
[G 0 -S r -G%(z,z,,z',£,£') 

PjjEz) 


(35) 


-h 


(E)V2ns l j(z — zi) 


exp • 


[E - (E".{z - z i ) )] 
2s"j{z-z\f 


Pk(E') 


<Pjk(( E 'k ( z i -z 7 ))) 




2[4(z-z') +E 2 jk \ 


|d£ 2 , 


(36) 


where ( E'!{z -z,)) = Rj 1 [/? y (£ 2 ) )] and 

Sj (z - zi) is the corresponding spread. The integral has a 
sharp maximum at E 2 = (E' k {z\ - z!))k jk = (E 2 ), E = 
( E'.(z-z \ )) = {Efiz)), where the cross-sections, attenua- 
tion functions, and straggling widths are evaluated. We 
expand {E"(z - zi)) about the maximal value of E 2 to 
obtain 


z(E"j (z - zi)) » (Efiz)) + r Jk [E 2 - ((E' k (z x - z')) - A^)], 

(37) 

where 

^ = [0^(z-z,)>] (Bi) =^Ml. (38 ) 

Substituting Eq. (37) into the integral (36), making 
the change of variables x = rj k [E 2 - ({E' k {zy - z / )) - /.jt)} 
and integrating with respect to x results in 


250 
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255 
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260 
261 
262 
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265 
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[C?-3'GP}j k (z,z u z’,E,E') 
P J ((E > t (z l -z'))-X jk )P k (E ') 
Pj{E)Pk{{E' k (z\ -/))) 

V / 2^) l ^ (2l) J 

271 where 


(39) 


fj(z 1 ) = 4 1 {*/[<4(z. - *0) - 4] - (z - z, ) } (40) 

273 and 

(*i) = y rjk [s*( z i - + 4] + s 'j( z ~ z < ) 2 - ( 4 1 ) 

275 Lastly, we need to evaluate the integral 

G) k (z,z\E,E') = [Q- G° • S r • G°] yi (z,z,,z' ,£,£') 

= J\(?-3'-(?} Jk (z,z u S,E,E')dz l . 

(42) 


277 For a given set of parameters z, £, z\ £', there is a value 

278 z„ of z\ at which the integrand of (42) achieves a max- 

279 imum and at which slowly varying factors entering the 

280 integrand may be computed. Thus 


282 

283 


Gjk(z,z',E,E') 


P J ((£' k (z m -z'))-A jk )P k (E') 
P J (E)P k ((E' i (z m -z'))) 
x ^(((ffUZm-zQ)) 
\J2rf k { Z m) 


X 



[E-fAzj)] 1 

240m) 


jdz,. (43) 


The point z m at which the integrand of (43) achieves its 
maximum is given by the equation 


fj{ z m) — E, 


(44) 


285 and is easily obtained by the routine root finding tech- 

286 niques. It is not difficult to show that 



5 t [(4(z w -z'))] 

Sj[(E' k {z m -z / ))-lj k ] 



(45) 


288 Therefore, in (43) we may use the substitution 

289 x = [E -/>(zi)]/[v / 2iy*(z m )] and then integrate to get 


C4(z,z\£,£') 


Pj[(E' k (z m -^)-X Jk }P k [E'\ 
P j {E)P k ((E' k (z m -z>))) 
a jk(((E' k {z m — Z 1 )))) 

X 2/;(z„) 

E-m 

V2s Jk {z m ) 


x < erf 


E-fjV) 

V2s Jk {z m ) 


- erf 


( 46 ) 


3.5. Second collision term 


292 


The third term in Eq. (31) is the second collision term 
G) k (z,zl,E,E') = [Q ■ G° • 2 r • G%{zJ,E,E') 

= [j G°.(z,z,,£,£ 2 l) 

x | Jlj °)p(Ei,E { ) 


X 


Gi t (zi,z',£’i,£ , ')d£i 


|d£ 2 dzi • 

J 


(47) 


293 


On making approximations similar to those used in the 295 
previous section, Eq. (47) is reduced to the form 296 



Fig. 1. Primary ion flux at various depths for Ne{20,10) incident on 
aluminum at 600 MeV/amu. 



291 where s jk (z, ) is given by (41). 


Fig. 2. First generation 0(16,8) fragment flux at various depths com- 
pared with previous results (broken curve). 
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G 2 (z z 1 E E') = V' f* - r (~E + )\ 

x G l pl! (z l ,z!,Ej + l Jp ,E')dz], (48) 

298 where 

£ y = /?7 1 [/? y (£)+z-z 1 ], (49) 

300 and is then evaluated by numerical quadrature. 


301 4, Results 

302 Shown in this section are some results for the gO 16 

303 fragments which are produced when a beam of ioNe 20 


ions strikes an aluminum target at 600 MeV/amu. The 304 
results presented are similar to those obtained for other 305 
ions. 306 

Fig. 1 shows the flux of the primary beam at various 307 
depths and exhibits the effects of energy straggling. In 308 
contrast to earlier works in which the primary beam 309 
appears as a propagating delta function, we see here that 3 1 0 

the primary beam attenuates and widens with depth. 311 

Note that the greatest depths in Fig. 1 are beyond the 312 
85% range where straggling propagators of the past 313 


have failed. 314 

Figs. 2 and 3 show the flux of the first and second 315 
generation of gO 16 ions produced, respectively, in com- 316 
parison with the corresponding flux (broken curves) 317 
obtained from an earlier approximate code using a non- 318 



Fig. 3. Second generation 0(16,8) fragment flux at various depths 
compared with previous results (broken curve). 



Fig. 5. Second generation 0(16,8) fragment flux at various depths 
compared with the same flux for the case in which the collision energy 
downshift is zero (broken curve). 




Fig. 4. First generation 0(16,8) fragment flux at various depths com- 
pared with the same flux for the case in which the collision energy 
downshift is zero (broken curve). 


Fig. 6. First generation 0(16,8) fragment flux at various depths com- 
pared with the same flux for the case in which the interaction energy 
width is zero (broken curve). 
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319 perturbative expansion of the solution (Wilson et al., 

320 1991). The non-perturbative approximation lacks spec- 

321 tral details and assumes a broad near uniform distri- 

322 bution over the allowed energy domain (Wilson et al., 

323 1991). 

324 The effect of the collision energy downshift is ex- 

325 hibited Figs. 4 and 5, where the flux of the first and 

326 second generation of g 0 5 * * * * * * * * * * 16 ions is compared with the 

327 corresponding results for the case in which the down- 

328 shift is zero. For the ions shown the shift is not great, 

329 contributing only a few MeV/nucleon. The downshift of 

330 more massive projectiles is somewhat larger. 

331 Figs. 6 and 7 exhibit the effect of the interaction en- 

332 ergy width £,* by comparing the flux of the first and 


Fig. 9. Second generation 0(16,8) fragment flux at various depths 
compared with the same flux for the case in which there is no energy 
straggling (broken curve). 


second generation of g O l6 ions with the corresponding 333 
results for the case in which the interaction energy width 334 
is zero. Significant widening occurs in the first genera- 335 
tion of secondaries but little effect is seen in the second 336 
and presumably higher generations. 337 

The effect of energy straggling on the the first and 338 
second generation of sO 16 ions is exhibited Figs. 8 and 9, 339 

where the flux of of these ions is compared with the 340 
corresponding results for the case in which no energy 341 
straggling is present. In contrast with the results for the 342 
primary beam where straggling makes a significant 343 
contribution the effect on the first and second generation 344 
of secondaries appears to be relatively small. 345 




E, MeV/rtudeon 



E, MeV/nudeon 


Fig. 8. First generation 0(16,8) fragment flux at various depths com- 
pared with the same flux for the case in which there is no energy 
straggling (broken curve). 


5. Concluding remarks 346 

The present formalism provides means of easy vali- 347 

dation of material transmission properties in conven- 348 

tional laboratory setups at least to the third 349 

perturbation term. Higher order terms can be easily 350 

added using non-perturbation theory and assuming the 351 

third term spectral distribution. The next step is the 352 

simulation of the detector responses so that “raw” ex- 353 

perimental data can be used to validate model predic- 354 

tions thereby simplyfying the validation process. Early 355 

versions of the Green’s function code, when coupled 356 

with multiple elastic scattering in terms of acceptance 357 
functions (Shavers et al., 1993), showed great promise in 358 
describing HZE ion transport. In those studies, the 359 
straggling, energy downshift, and dispersion were ne- 360 
glected. The present formalism corrects those last re- 361 
maining deficiencies. The recognition of the present 362 
formalism as the lowest order asymptotic term provides 363 
a systematic approach to more realistically treat a host 364 
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of ion beam related problems. The next step will be to 
couple the multiple scattering propagator to the for- 
malism and adding transverse momentum components 
to the first interaction term. 


6. Uncited references 

Schimmerling et al. (1999); Tschalar and Maccabee 
(1968). 
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